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ABSTRACT 

The Sersic model is the de facto standard to describe the surface brightness distribution of hot stellar systems. An important incon- 
venience of this analytical model is that the corresponding luminosity density and associated properties cannot be expressed using 
elementary functions or even standard special functions. We present a set of compact and elegant analytical expressions for the lu- 
minosity density, cumulative luminosity and potential for the Sersic model in terms of the Fox H function for general values of the 
Sersic index. Furthermore, we present explicit series expansions of these quantities and discuss the asymptotic behaviour. Our analysis 
completes the work of Mazure & Capelato (2002) and Baes & Gentile (201 1) and demonstrates the power of the underestimated Fox 
H function as a tool for analytical work. 

Key words, methods: analytical - galaxies: photometry 



1. Introduction 

The Sersic (1968) surface brightness profile, has become the preferred model to describe the surface brightness profile of early- 
type galaxies and the bulges of spiral galaxies (e.g. Davies et al. 1988; Caon et al. 1993; D'Onofrio et al. 1994; Cellone et al. 1994; 
Andredakis et al. 1995; Prugniel & Simien 1997; Mollenhoff & Heidt 2001; Graham & Guzman 2003; Allen et al. 2006; Gadotti 
2009). Many analytical properties of this model have been discussed in the literature (e.g. Ciotti 1991; Ciotti & Lanzoni 1997; 
Ciotti & Bertin 1999; Trujillo et al. 2001; Mazure & Capelato 2002; Cardone 2004; Graham & Driver 2005; EKasdottir & Moller 
2007; Baes & Gentile 2011). 

An important inconvenience of the Sersic model is that its deprojected luminosity density, i.e. the spatial 3D luminosity density 
v(r) that projects on the plane of the sky to the Sersic surface brightness profile, cannot be expressed using elementary functions or 
even in terms of standard special functions. It was long thought that no analytical expression could be obtained. Quite unexpectedly, 
Mazure & Capelato (2002) came up with an analytical expression for v(r) in terms of the Meijer G function for all integer Sersic 
indices m. Baes & Gentile (201 1, hereafter BG1 1) took this analysis one step further and showed that the deprojection of the Sersic 
surface brightness profile for general values of m can be solved elegantly using Mellin integral transforms and gives rise to a Mellin- 
Barnes integral. The result is that the Sersic luminosity density can be written compactly in terms of a Fox H function, which 
reduces to a Meijer G function for all rational values of m. BG11 also calculated analytically a number of additional properties of 
the Sersic model for rational m, including the asymptotic expansion of the luminosity density at small and large radii, the cumulative 
light profile and the gravitational potential. This was possible thanks to the many analytical properties of the Meijer G function. 

The goal of this paper is two-fold. Foremost, we extend and complete the analysis presented in Mazure & Capelato (2002) and 
BG11: we will provide compact and elegant expressions for the density, potential and luminosity profiles in terms of the general 
Fox H function, which are valid for all values of the Sersic index m rather than just for integer or rational m. We also present a 
completely general series expansion of these functions that enables both a numerical evaluation and a straightforward analytical 
study of the asymptotic behaviour. Besides providing these useful characteristics on one of the most used models in extragalactic 
astronomy, our work also has a secondary objective, namely to demonstrate the power of the Fox H function for analytical work. 
The Fox H function is, to use an understatement, not the most mainstream special function: it is even not contained in the standard 
works on special functions such as Gradshteyn & Ryzhik (1965) or the Wolfram Functions Site. We feel that this is not justified: it 
is in fact a very elegant and powerful tool for analytical work and it is becoming more and more used in mathematics and applied 
sciences, including physics, biology, engineering and earth sciences. It is also gradually being used in astrophysics, e.g. for solar and 
stellar structure models, fractional reaction-diffusion equations and stellar dynamics (Haubold et al. 2007, 2011; Baes & van Hese 
2007; Van Hese et al. 2009; De Rijcke et al. 2010). With our analysis, we wish to illustrate its useful properties and advocate its use 
in theoretical astrophysical research. 

In Section 2 we derive compact expressions for the luminosity density, cumulative luminosity and potential for the Sersic model, 
and in Section 3 we deduce detailed power and logarithmic-power series expressions for these important quantities. We discuss the 
asymptotic behaviour of our functions in Section 4, and we sum up our results in Section 5. In Appendix A we present the Fox H 
function and we discuss some of its properties that were used for our analysis. 
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2. Analytical properties of the Sersic model 

The Sersic model is defined by the intensity profile projected on the plane of the sky, 
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The 3D, deprojected luminosity density v(r) of a spherically symmetric system can be recovered from the surface brightness profile 
I(R) using the standard deprojection formula 
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Substituting the Sersic profile (1) into (2) we obtain an integral that cannot readily be evaluated using the standard ways or look-up 
tables. BG1 1 applied a Mellin integral transform technique to convert this integral to a Mellin-Barnes contour integral, 
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where we have used the reduced coordinate 
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BG11 converted this general expression (or actually, a slightly different but equivalent expression) to expressions in terms of the 
Meijer G function for integer and rational values of the Sersic index m. As a check on these formulae, they calculated the total 
luminosity of the Sersic model for rational m by integrating the luminosity density over the entire space. To obtain the results, they 
used the integration properties of the Meijer G function, combined with several applications of Gauss' multiplication theorem. As 
a generalization of this result, and as a nice example of the power of the Fox H function, we calculate the total luminosity from the 
general formula (4), i.e. 
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To evaluate this integral, one should remind that equation (A.l) defined the Fox H function as an inverse Mellin transform of a 
combination of gamma functions. As a result, the Mellin transform of a Fox H function reads 
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Applying this to (6) with s - 1, we obtain 
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in agreement with the value obtained by integrating the surface brightness profile (1) over the plane of the sky (Ciotti 1991). 

From the luminosity density, a number of other important quantities can be derived, most importantly the cumulative luminosity 
profile L(r) and the gravitational potential ^(r), 
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where the Y is the mass-to-light ratio. Mazure & Capelato (2002) and BG11 calculated these quantities for the Sersic model for 
integer and rational values of the Sersic parameter, respectively, using the integration properties of the Meijer G function. It is, 
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however, possible to calculate these properties for general m in an elegant way by directly applying the integrations on the Mellin- 
Barnes integral form of the luminosity density. For the cumulative luminosity profile we find 
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For the gravitational potential we find after a similar calculation 
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The formulae (4), (11) and (12) form a triplet of formulae that describe three important spatial properties of the Sersic model in a 
compact way. 

A straightforward way of checking these formulae is to look at the model that corresponds to m = |. In this case, all components 
of the vectors A and B are equal to one, such that the Fox H function reduces to a Meijer G function. We find 
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These expressions can also be derived by substituting the intensity profile I(R) - lot bR 'l R l in the expressions (2), (9) and (10) and 
directly evaluating the resulting integrals. More generally, one can check that the formulae (4), (11) and (12) reduce to the equations 
(22), (40) and (44) of BG1 1 for rational values of m. 



3. Explicit series expansions 

While the expressions (4), (11) and (12) form an triplet of compact formulae that are useful for analytical work, they are not readily 
useful to numerically evaluate the spatial properties of the Sersic model. For rational values of m, the Fox H functions reduce to 
Meijer G functions, and some numerical software packages have this function now implemented. However, the numerical evaluation 
of Meijer G functions with large parameter vectors (which easily occurs in our case for rational values of m, as can be seen in BG1 1 ), 
proves to be difficult, in particular in cases where second-order poles are present in the integrand of the inverse Mellin transform. 
Moreover, for general values of m, the expressions (4), (11) and (12) cannot be written in terms of the Meijer G function or any other 
special function, and we are not aware of any implementations in numerical software that can evaluate general Fox H functions. 

In this section, we derive explicit series expansions for v(r), L(r) and ^(r), which both enables the numerical calculation and 
again highlights the power of Fox H function as a useful mathematical tool. The expansions build on the general series expansion 
of the Fox H function as power or power-logarithmic series (details can be found in Appendix A). The form of the series expansion 
depends on the multiplicity of the poles of the gamma functions T(bj + Bjs). For v(r) and L(r), the poles of these gamma functions 
are found at -k[/2m and -1/2 — hi with k\ and ki any natural number. The gamma functions corresponding to the expression of 
the potential T(r) contain the same poles with an additional pole at 1 /2. The good news is that each pole can at most occur twice, 
the bad news is that this happens quite often: for all integer m and rational m = p/q where the denominator q of the fraction is odd, 
double poles do occur. 

In the case of simple poles of the gamma functions Y(bj + Bjs), the expansion of the Fox H function is a power series, given by 
equation (A. 6). Applied to our case, we find in case m is non-rational or rational with an even denominator, the expansions 
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Notice that the term k — in the first sums is omitted, since the factors F(0) in the denominator make those terms vanish. In fact, if 
m is a rational number p/q with q even, then the terms in the first sums for which k — 0, p, 2p, . . . vanish; if p = 1, these first sums 
vanish completely. A particularly interesting case is (again) m=\, where we find 
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in agreement with formulae (13), (14) and (15). 

When the gamma functions T(bj+Bjs) have multiple poles, the expansion of the Fox H function is a logarithmic-power series, the 
complexity of which increases with increasing multiplicity of the poles. The full expression for the case where two gamma functions 
share poles is given in equation (A. 10). For our present case, this means that we obtain logarithmic -power series expansions for 
v(r), L(r) and *F(r) if the Sersic index m is integer or rational with an odd denominator. If we define k$ — (q + I) 12, one obtains after 
quite some algebra 
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In these expressions, if/(s) = r'(s)/r(s) is the digamma function. In ^(r), the third term 2T(m)ly[n corresponds with the residue of 
the pole I. Again, the terms in the first sums for which k = 0,p,2p, . . . vanish, and these first sums vanish completely if p = 1, 
On the other hand, if q = 1, i.e. if m is an integer value, the second sums vanish, since then the integrands for v(r) and L(r) 
have no simple poles /?2,,t, , while for ¥(r) only the pole I remains as a single pole. On the crossroad of these two cases we have 
m = 1: for this model, both the first and second sums in the expansions (22)-(24) vanish completely, apart from a single term for 
the potential. Astrophysically, the Sersic model with m = 1 corresponds to a model with an exponential surface brightness profile, 
often used for the description of low-luminosity elliptical galaxies and pseudo-bulges. BG11 calculated the luminosity density of 
the exponential model by directly deprojecting the surface brightness profile and through its representation as a Meijer G function. 
For the luminosity density we get 
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with K v (z) the modified Bessel function of the second kind. This expression is in agreement with equation (24) of BG1 1. Similarly, 
we obtain for the cumulative luminosity 



tf- z-j (2k + 3) k\ W 



-In 

2tt/ ^ 

b 2 



(-)+ —^—+il/(k+\) 
\2J 2k + 3 V 



Ki(u)L\ (u) + K\(u) Liiu) uK\(u) 

3n 



with L v (u) the modified Struve function. Finally, for the potential of the exponential model we get the expansion 
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These last two expressions can also be obtained by substituting the luminosity density (25) into the recipes (9) and (10). 
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4. Asymptotic behaviour 

With all the explicit power series of Section 3 available, it is fairly straightforward to examine the asymptotic behaviour of the 
spatial function of the Sersic model at small radii, generalizing the results of Ciotti (1991) and BG11. The density has the following 
rich behaviour, depending on the value of m: 
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with y a 0.57721566 the Euler-Mascheroni constant. The luminosity behaves as 
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Finally, the potential approaches r — > as 
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If we set r = in the expressions (30), we recover the central potential 
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for every real value m > 0, in agreement with equation (12) of Ciotti (1991). 
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5. Conclusions 

This paper is a companion paper to BG11 that deals with the derivation of analytical expressions for the luminosity density 
v(r), cumulative luminosity L(r) and gravitational potential *F(r) of the Sersic model. Our work extends the work initiated by 
Mazure & Capelato (2002), who managed to express these properties in terms of the Meijer G function for integer values of the 
Sersic index m. In BG11 we extended this analysis by demonstrating that these expressions can be derived using a Mellin integral 
transform and that they are also valid for half-integer values of m. We also derived more general expressions in terms of the Meijer 
G function for these spatial properties for Sersic models with rational m. Actually, the Mellin integral transform approach directly 
led to an expression in terms of the Fox H function for the luminosity density valid for general values of m. Not fully aware of the 
rich power of the Fox H function and attracted by the availability of implementations of the Meijer G function both in symbolic 
computer algebra packages and as high-performance computing code, we focused our attention in BG1 1 on the study of the Sersic 
model with integer and rational m. In this paper, we have extended this analysis to arbitrary values of m, now fully making use 
of the Fox H function. We have derived compact and elegant expressions for the luminosity density, cumulative luminosity and 
gravitational potential that are valid for all values of m. We have used the properties of the Fox H function to deduce explicit power 
and logarithmic-power expansions for v, L and *P that can be used to evaluate these important spatial properties for arbitrary values 
of m. These series expansions also provide a direct way to probe the rich asymptotic behaviour of the Sersic model at small radii; 
our results are in full agreement with the more ad hoc analysis of BG1 1. 

Our work fits into a significant effort of analytical work on the Sersic model: apart from the spatial properties discussed 
here, also its photometric (Graham & Driver 2005), dynamic (Ciotti 1991; Ciotti & Lanzoni 1997) and lensing (Cardone 2004; 
EKasdottir & Moller 2007) properties have studied extensively through analytical means. Given that the Sersic model is the de facto 
standard model to describe the surface brightness distribution of hot stellar systems, this combined range of analytical properties has 
plenty of applications. For example, they can be used for the construction of self-consisting dynamical models for elliptical galaxies, 
serve as a realistic starting point for numerical A^-body or smoothed particle hydrodynamics simulations, or serve as a template for 
radiative transfer or gravitational lens modelling. We are convinced that the exact analytical results for the spatial properties of the 
Sersic models will be a useful contribution in this respect. 

Besides this practical use, we hope to have demonstrated the power of the Fox H as a tool for analytical work. We are aware 
that the Fox H function is not the most every day's special function and that its unfamiliar definition through inverse Mellin 
transform alone might already seem daunting. However, the Fox H function is gradually appearing more in mathematics and applied 
sciences and there are now several volumes of comprehensive literature available on its properties (e.g. Mathai & Saxena 1978; 
Srivastava et al. 1982; Kilbas & Saigo 2004; Mathai et al. 2009). Through this work, we advocate the use of the Fox H function in 
theoretical astrophysics. 
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Appendix A: The Fox H function 



Tjin.n 
H P.q 



(a, A) 




(b,B) 


z 



The Fox H function, or just called the H function, is generally defined as the inverse Mellin transform of a product of gamma 
functions, 

1 r Il%iT(bj + B i sm%im-aj-Ajs) 

— d z ds - ( A -i) 

2m Jl U% m+1 T(l - bj - Bjs) U%„ +i r(aj + Ajs) 

It was introduced by Fox (1961) in his attempt to find the most general symmetric Fourier kernel. The Fox H function is a gen- 
eralization of the Meijer G function, but includes also many other special functions, including Mittag-Leffler functions and gener- 
alized Bessel functions. For details on the definition, convergence and many useful properties of the Fox H function, we refer to 
Mathai & Saxena (1978), Srivastava et al. (1982), Kilbas & Saigo (2004) or Mathai et al. (2009), and the references therein. 
We will now derive a series expansion of the Fox H function. Let us write the Fox H function in the form 
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with 



<p(s) = 



= ^ J <f(s)z s ds, 
U7 = inbj + Bjs)U n j=i ni-aj-Aj S ) 

U% m+l m - bj - B jS ) n% n+1 r(aj + A jS ) • 



(A.2) 



(A.3) 



Under certain conditions (which are always satisfied for the Fox H functions we consider in this paper; see e.g. Kilbas & Saigo 
(1999) or Mathai et al. (2009) for details), the Fox H function is an analytical function, and the contour integral can be evaluated 
using the residue theorem. If we introduce for the poles of the functions T(bi + Bjs) the short-hand notation = -(b, + k)/B,, then 



Hi 



(a, A) 




(b,B) 


z 



(=1 k=0 

If B, tji is a simple pole, then the corresponding residue is fairly straightforward: we find 



r _ I _ (-1)* Z~ Pl * 

Res Ms) z s \ = lim (s - Bt,k) f(s) z " = tfii(B uk ) — — , 



(A.4) 



(A.5) 



s=Pu 1 ' ' ' J ' " "" """ k\ Bj 

where (fi(s) is ip(s)/T(bi + BfS). As a result, if all gamma functions T(bi + B,s) have only single poles, we get as a series expansion, 
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Jbi+k)/Bi 



(A.6) 



in agreement with equation (3.12) in Kilbas & Saigo (1999). However, if several gamma functions share the same pole, then this 
pole is of second order, and the calculation becomes more cumbersome. Kilbas & Saigo (1999) demonstrate that the Fox H function 
can then be expressed as a logarithmic-power series rather than a simple power series. They present a generic expression valid for all 
orders of pole multiplicity. Here, we present a less general, but more explicit, series expansion in the case that two gamma functions 
Y(bj + Bjs) share at least one pole. Without loss of generality, we can place these two gamma functions at the front, so that we can 
write 

<f(s) = r(bi + B x s) T(b 2 + B 2 s) <f>( S ), (A.7) 

with 

U7 =3 nbj + Bjsm n j =1 T(l-aj-Ajs) 



<t>(s) 



n% m+l ni-bj-B jS )n p j =n+1 T(aj + A jS y 



(A.8) 
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Now, suppose there is a tuple (k\ , k 2 ) of indices so that = B 2 jc 2 ■ The residue of this second-order pole is then, after some algebra, 



Res fo(*)z""'} = lim (1 Us -B lM f <p(s) Z - s ] 



•lnz + Bi^i + l) + B 2 if/(k 2 + 1) 
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(A.9) 



with i//(s) = F'(s)/r(s) the digamma function. Moreover, note that (p(s) is a product and quotient of gamma functions, so that 
(p'(s)/(p(s) can also be expressed as a sum of digamma functions. This means that all the machinery is available to express the Fox 
H function as the rather daunting-looking series expansion 
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with the constants Q. defined as 



Jbi+kO/B, 



(Q.-lnz), (A. 10a) 



C kl = Bnj/(k\ + 1) + B 2( A(A:2 + D + J] ^'-A - Bffi) -J]A^(l - a } + Aj b -^) 

7=1 

+ J B^(l-b j + Bj^)- J] A^-A,^). (A.lOb) 



7=3 



/=m+l 



j=n+l 



The prime in the first summation in (A. 10a) indicates that this sum covers only the single poles, and the double prime in the second 
summation indicates that this summation runs over the second-order poles. In the latter summation, we set k 2 = B 2 (b\ +k l )/B 1 - b 2 . 



